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Abstract 

The statistics of outbreaks in a model for the propagation of meningococcal diseases 
is analyzed, taking into account the possibility that the population is fragmented 
into weakly connected patches. It is shown that, depending on the size of the sample 
studied, the ration between the variance and the average of infected cases can vary 
from unity (Poisson statistics) to e~^, where e is the normalized infection rate. 
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1 Introduction. 



The meningococcus is a major cause of meningitis and septicaemia. Despite 
this, infection with the meningococcus is mostly harmless and only rarely leads 
to disease. Transmission of the disease is almost exclusively through asymp- 
tomatic carriers of the disease. A predominant feature of the epidemiology 
of meningococcal disease are outbreaks of variable scale and duration. The 
meningococcal population is genetically highly diverse. We have shown, us- 
ing a mathematical model, that heritable diversity with respect to pathogenic 
potential can lead to disease outbreaks [1,2,3]. 

Meningococcal disease is a a notifiable disease in many countries. Therefore 
there exist extensive data sets on the incidence of meningococcal disease. The 
analysis of meningococcal disease data is problematic because the number of 
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asymptomatic carriers at any time, the variable that is probably of most inter- 
est, is normally not known because transmission of the pathogen takes place 
almost exclusively through asymptomatic carriers. Therefore key epidemio- 
logical parameters are difficult to estimate and methods that are standard in 
epidemiology, such as outbreak reconstruction through contact tracing, can 
not easily be applied. For this reason outbreaks of meningococcal disease are 
difficult to reconstruct and to detect. 

In this paper we will investigate the statistical structure of an epidemiological 
model to infer the underlying disease process from data on the number of cases 
of disease. Such insights have been applied in the analysis of meningococcal 
disease data [3]. Here we will investigate the validity of the assumptions made 
for these inferences and study how the variance in the number of cases of 
disease depends on the structure of the population. 



2 The SIRYX model. 



We study the SIRYX model, considered in [1,2,3]. The model is an extension of 
the SIR model[4]. There are two types of infected individuals, / and Y. The Vs 
are generated by mutation from the J's at rate fif3. For simplicity we assume 
that the back mutation rate Y ^ I is nil. The Y population can develop 
disease at rate e/3. The parameter e is the pathogenicity: the probability to 
develop disease upon infection. We define the number of individuals which 
suffer the disease X. We further simplify the model by assuming that these 
individuals are removed from the population. The mean field equations are: 



The only difference with respect to the model studied in[l,3] is the intro- 
duction of a the rate 6 at which the X individuals are removed from the 
population (see below). This rate implies, that, in the long run, the only sta- 
tionary situation is the conversion of all individuals into the X type, and their 
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eventual disappearance. We will study here quasisatationary situations which 
arise when S,e <^ 1. 

Following [1], we consider that the system is near its stable point when /i — 
and Y — 0. The remaining parameters at the fixed point are: 



S__S^_ 7 



Assuming that the fixed point values for S and / do not change much for small 

/X, wc can define a simple birth-death model for the variables Y and X. We 
define p{Y^t) = J2x=oPO^' ^ probability of finding the value Y at 

time t. This function satisfies: 



-p{Y,t)^[b{Y-l) + c]p{Y-l,t) 

+ a{Y + l)p{Y + l,t)- {bY + aY + c)p{Y, t) (3) 

where we have defined the death rate, a = 7, birth rate, b = (3{1 — e)S*/Np 
and rate of creation of a new individuals by mutation, c = Pii^S* /Np)I* . 

We generahze this equation to the case of a system divided into M patches. 
The main difference is that the birth probability has to be divided into the 
probability that the contagion is to another individual within the same patch, 
which we still define as 6, and the probability that the contagion leads to a 
new individual of type Y in another patch, b' . The total infection rate remains 
6 + (M - 1)6' = - e)S*/Np. The generahzation of eq. (3) is: 



-p{{Y,} , = ^ b{Y, - l)p{Y, - 1, {Yj} , + E E b%PiY^ - 1, {Y^ , t) + 

+ E «(^^ + my^ + 1, {Yj} , ^) + E - 1' {^i} ' - 

i i 

- E bY,p{{Y^ , ^) - E E b'Yjp{{Y^t , t) - 
-Y,aYip{{Yi},t)-cp{{Y,},t) (4) 

i 

Note that in this equation c is the mutation rate within one patch. The to- 
tal mutation rate is Mc. When b = b' we recover the limit of a well mixed 
population, while for b' — the patches are decoupled. 
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3 Results. 



Prom eq.(4) we can calculate the ensemble means of different quantities. The 
details of the calculations are given in the Appendix. The results are: 



{Yd 



a-b-{M - l)b' 

{y^)-{yi) 



2 2 «c[a -b-{M - 2)6'] 



[a-b- {M -l)bf{a-b + b') 
(YiVj) - {Yi}{Yj) = (5) 

All these quantities vanish when the mutation rate is zero, c = 0. 

The net growth rate is ae ~ a — b — (M — 1)6'. The number of infected cases 
appear with rate X'- — aeYj,. 

We now calculate the number of infected individuals, X^. We study first the 
case of a single population and a single variable X. The infected individuals 
are generated from the Y^s at rate aa — b — (M — 1)6'. In order to calculate X 
in a single population, we use as unit of time a^^, and assume that the death 
rate of the X's is 6. We write the mutation rate c = aj. then, we can write: 



^^^^^^M ^e[(Y + 1)P{Y +l,X-l,t)- YP(X, Y, t)] + 

+ 5[{X + 1)P{X + 1, Y, t) - XP{X, Y, t)] - jP{X, Y, t) (6) 

where we have used as the unit of time a~^, e is now the rate of conversion 
from Y into X, 5 is the death rate of the X's, and 7 = c/a is the mutation 
rate from I into Y. Using the techniques described in the Appendix, we can 
write: 



dt 



d{X 
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: 2e{XY) + e{Y) - 25{X') + 5{X) 



dt 

^ -e{XY) + e{Y'') - e{Y) - 6{XY) + c{X) (7) 



d{XY) _ , ^,^^2 



dt 

In a stationary state the right hand side of these equations is equal to zero, 
and we find: 



{X)^'^{Y) 
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e + 

We substitute the first and third of these equations into the second, so that: 



Prom eq.(6) we also obtain (Y) = 7/e. Inserting this result into eq.(9), we 
have: 

(X) 2 e + 5 e + 5 (Y) ^ ' 



This equation relates the variance and the average of X. When the mortality 
rate is very high, 5 ^ e, we have: 



(11) 



5>e 



The ratio approaches a constant of order unity, and the process seems to 
have Poisson statistics. This is reasonable, because there is an approximately 
constant reservoir of Y individuals which can lead to an X individual which 
disappears quickly, and the distribution of X cases is not influenced by the 
fluctuations of Y. 

A more interesting regime arises if 5 <^ e and e <^ 1. Then, the r.h.s. in eq.(lO) 
is dominated by the third term, because ((F^) — ~ We find in 

this case: 



This result is the basis of the following section. Note that when 5 = the 
value of {Xi) increases linearly with time. 
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Fig. 1. Dependence of the variance of infected individuals on the average for small 
populations which are part of a larger, well mixed population, {X) < {X*), or form 
an isolated population, (X) > (X*) (see eq.(17)). The infection rate is e = 0.01. 

4 Size effects. 



Using the results in the Appendix and eq.(12) we find (for M > 2): 

VaiXj _ a[a-b-{M- 2)b'] 
{Xi) ~ [a-h-{M -l)h']{a-h + V) 



(13) 



On the other hand, for the entire system we obtain: 
VarX a 



{X) a-b-{M - l)b' 



(14) 



The linear relationship between the variance and the mean is discussed in 
detail in[3]. For isolated patches, b' = and 6/a = 1 — e. As expected, the 
local and global values, eq.(13) and eq.(14) coincide, giving a ratio equal to 

1/6. 

In a well mixed population, we have b' = b, the total birth rate is 6tot — 

Mb, 

and (Mb) /a = 1 — e. Then, we find: 



VarX^ 1-e 

VarX 1 , , 
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For a small subsystem of a well mixed population (M 1/e), wc have 
VarXj/(Xj) ^ 1. This ratio would imply that the process is due to random 
mutations with Poisson statistics. An analysis of the total variance, however, 
gives a rather different result. For large (but artificial) subdivisions of the well 
mixed system. Me < 1, and YoxXi/ {Xi) ^ l/(Me). 

It is interesting to analyze the situation in which populations of size N below 
some size N* are part of a well mixed population of size iV*, while larger 
populations can be considered as isolated, made up of smaller, decoupled pop- 
ulations of size A^*. Then, for populations < A"* we can use eq.(15) with 
M = N*/N = {X)/{X*) ((X*) is the value of the mean of a population of 
size A^*), while when A^ > A^* we can use eq.(16). The variance can be written 
as: 



. Eq.(17) interpolates between a Poisson like regime for A^ -C A^* to a 1/e 
ratio between the variance and the mean for A" > A"*. A sketch of the results 
is shown in Fig.[l]. used here imply that the coupling between different parts 
of the 
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(17) 



^{b-a){Y^ + b'Y.{yj) + c 



dt 



= 2(6 - a) ((r/) - + {a + bm + c + 2b' j:{Yj) + b' ^{Yi) 



dt 



d{Y^,)-{Yi){Y,) 



= 2(6-a)((F,F,)-(F,)(F,)) + 

+ b'Y: {Yiyk) - (Y^iYk) + (YjYk) - {Yj){Yk) + 



dt 
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so that limt^oo{Yi) = Ke-^''-''-^'^-^^'''^^ + c/[a - b - {M - 1)6'] where X is a 
constant determined by the initial conditions. We define: 



C,,^{Y,Y,)-{Y,){Y^) 
These quantities satisfy: 



d I Cii 



+ 



^2(6 -a) 2(M-1)6' 
^ 2h' 2{b - a) + {M - 2)b' ^ 

a + b){Y^+b'EjMYj) + c 




Ci', 



+ 



(.3) 



At long times, we find: 



lim 

t—i-OO 



{a + b){Y^it) + b'J2{Yj){t) + c 



2ac 



a-b-{M- 1)6' 



(-4) 



This leads to: 



lim 



2(6 -a) 2(M-1)6' 
26' 2(6 - a) + (M - 2)6' 



-1 


/ 2ac \ 


1 


a-b-(M-l)b' 




[ J 



:-5) 
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